{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# B 14-2 生存率曲线及比较\n",
    "\n",
    "Kaplan-Meier 生存率曲线绘制及 log-rank 检验。\n",
    "\n",
    "## 案例\n",
    "\n",
    "研究者收集了 15 例乳腺恶性肿瘤直径 <= 2cm 的患者和 21 例肿瘤直径大于 5cm 患者手术后的生存资料，定义生存时间为手术后到患者复发的时间为生存时间（月）。\n",
    "\n",
    "### 数据"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div><style>\n",
       ".dataframe > thead > tr,\n",
       ".dataframe > tbody > tr {\n",
       "  text-align: right;\n",
       "  white-space: pre-wrap;\n",
       "}\n",
       "</style>\n",
       "<small>shape: (36, 6)</small><table border=\"1\" class=\"dataframe\"><thead><tr><th>group</th><th>index</th><th>t</th><th>d</th><th>c</th><th>n</th></tr><tr><td>str</td><td>i64</td><td>i64</td><td>i64</td><td>i64</td><td>i64</td></tr></thead><tbody><tr><td>&quot;2cm&quot;</td><td>1</td><td>10</td><td>1</td><td>0</td><td>15</td></tr><tr><td>&quot;2cm&quot;</td><td>2</td><td>10</td><td>0</td><td>1</td><td>14</td></tr><tr><td>&quot;2cm&quot;</td><td>3</td><td>13</td><td>1</td><td>0</td><td>13</td></tr><tr><td>&quot;2cm&quot;</td><td>4</td><td>18</td><td>1</td><td>0</td><td>12</td></tr><tr><td>&quot;2cm&quot;</td><td>5</td><td>25</td><td>0</td><td>1</td><td>11</td></tr><tr><td>&hellip;</td><td>&hellip;</td><td>&hellip;</td><td>&hellip;</td><td>&hellip;</td><td>&hellip;</td></tr><tr><td>&quot;5cm&quot;</td><td>17</td><td>47</td><td>1</td><td>0</td><td>5</td></tr><tr><td>&quot;5cm&quot;</td><td>18</td><td>52</td><td>1</td><td>0</td><td>4</td></tr><tr><td>&quot;5cm&quot;</td><td>19</td><td>54</td><td>1</td><td>0</td><td>3</td></tr><tr><td>&quot;5cm&quot;</td><td>20</td><td>60</td><td>1</td><td>0</td><td>2</td></tr><tr><td>&quot;5cm&quot;</td><td>21</td><td>86</td><td>1</td><td>0</td><td>1</td></tr></tbody></table></div>"
      ],
      "text/plain": [
       "shape: (36, 6)\n",
       "┌───────┬───────┬─────┬─────┬─────┬─────┐\n",
       "│ group ┆ index ┆ t   ┆ d   ┆ c   ┆ n   │\n",
       "│ ---   ┆ ---   ┆ --- ┆ --- ┆ --- ┆ --- │\n",
       "│ str   ┆ i64   ┆ i64 ┆ i64 ┆ i64 ┆ i64 │\n",
       "╞═══════╪═══════╪═════╪═════╪═════╪═════╡\n",
       "│ 2cm   ┆ 1     ┆ 10  ┆ 1   ┆ 0   ┆ 15  │\n",
       "│ 2cm   ┆ 2     ┆ 10  ┆ 0   ┆ 1   ┆ 14  │\n",
       "│ 2cm   ┆ 3     ┆ 13  ┆ 1   ┆ 0   ┆ 13  │\n",
       "│ 2cm   ┆ 4     ┆ 18  ┆ 1   ┆ 0   ┆ 12  │\n",
       "│ 2cm   ┆ 5     ┆ 25  ┆ 0   ┆ 1   ┆ 11  │\n",
       "│ …     ┆ …     ┆ …   ┆ …   ┆ …   ┆ …   │\n",
       "│ 5cm   ┆ 17    ┆ 47  ┆ 1   ┆ 0   ┆ 5   │\n",
       "│ 5cm   ┆ 18    ┆ 52  ┆ 1   ┆ 0   ┆ 4   │\n",
       "│ 5cm   ┆ 19    ┆ 54  ┆ 1   ┆ 0   ┆ 3   │\n",
       "│ 5cm   ┆ 20    ┆ 60  ┆ 1   ┆ 0   ┆ 2   │\n",
       "│ 5cm   ┆ 21    ┆ 86  ┆ 1   ┆ 0   ┆ 1   │\n",
       "└───────┴───────┴─────┴─────┴─────┴─────┘"
      ]
     },
     "execution_count": 1,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "import polars as pl\n",
    "\n",
    "df_raw = pl.read_csv(\"B_14_2-data.csv\")\n",
    "\n",
    "df = df_raw.clone()\n",
    "\n",
    "df"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 生存分析"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "shape: (15,)\n",
      "Series: 't' [i64]\n",
      "[\n",
      "\t10\n",
      "\t10\n",
      "\t13\n",
      "\t18\n",
      "\t25\n",
      "\t…\n",
      "\t54\n",
      "\t68\n",
      "\t71\n",
      "\t88\n",
      "\t95\n",
      "]\n"
     ]
    }
   ],
   "source": [
    "from lifelines import KaplanMeierFitter\n",
    "\n",
    "T2 = df.filter(pl.col(\"group\") == \"2cm\").get_column('t')\n",
    "E2 = df.filter(pl.col(\"group\") == \"2cm\").get_column('d')\n",
    "T5 = df.filter(pl.col(\"group\") == \"5cm\").get_column('t')\n",
    "E5 = df.filter(pl.col(\"group\") == \"5cm\").get_column('d')\n",
    "\n",
    "kmf2 = KaplanMeierFitter()\n",
    "kmf5 = KaplanMeierFitter()\n",
    "\n",
    "print(T2)\n",
    "\n",
    "kmf2.fit(T2, E2)\n",
    "kmf5.fit(T5, E5);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 生存率\n",
    "\n",
    "乳腺肿瘤直径 $ \\leqslant $ 2cm :"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>removed</th>\n",
       "      <th>observed</th>\n",
       "      <th>censored</th>\n",
       "      <th>entrance</th>\n",
       "      <th>at_risk</th>\n",
       "      <th>KM_estimate</th>\n",
       "      <th>KM_estimate_lower_0.95</th>\n",
       "      <th>KM_estimate_upper_0.95</th>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>event_at</th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0.0</th>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>15</td>\n",
       "      <td>15</td>\n",
       "      <td>1.000000</td>\n",
       "      <td>1.000000</td>\n",
       "      <td>1.000000</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>10.0</th>\n",
       "      <td>2</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>15</td>\n",
       "      <td>0.933333</td>\n",
       "      <td>0.612641</td>\n",
       "      <td>0.990332</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>13.0</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>13</td>\n",
       "      <td>0.861538</td>\n",
       "      <td>0.549702</td>\n",
       "      <td>0.963561</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>18.0</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>12</td>\n",
       "      <td>0.789744</td>\n",
       "      <td>0.479128</td>\n",
       "      <td>0.927070</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>25.0</th>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>11</td>\n",
       "      <td>0.789744</td>\n",
       "      <td>0.479128</td>\n",
       "      <td>0.927070</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>29.0</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>10</td>\n",
       "      <td>0.710769</td>\n",
       "      <td>0.398260</td>\n",
       "      <td>0.881081</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>30.0</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>9</td>\n",
       "      <td>0.631795</td>\n",
       "      <td>0.325863</td>\n",
       "      <td>0.828573</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>33.0</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>8</td>\n",
       "      <td>0.552821</td>\n",
       "      <td>0.260079</td>\n",
       "      <td>0.770388</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>46.0</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>7</td>\n",
       "      <td>0.473846</td>\n",
       "      <td>0.200256</td>\n",
       "      <td>0.706899</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>50.0</th>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>6</td>\n",
       "      <td>0.473846</td>\n",
       "      <td>0.200256</td>\n",
       "      <td>0.706899</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>54.0</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>5</td>\n",
       "      <td>0.379077</td>\n",
       "      <td>0.129991</td>\n",
       "      <td>0.630543</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>68.0</th>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>4</td>\n",
       "      <td>0.379077</td>\n",
       "      <td>0.129991</td>\n",
       "      <td>0.630543</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>71.0</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>3</td>\n",
       "      <td>0.252718</td>\n",
       "      <td>0.049287</td>\n",
       "      <td>0.533373</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>88.0</th>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>2</td>\n",
       "      <td>0.252718</td>\n",
       "      <td>0.049287</td>\n",
       "      <td>0.533373</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>95.0</th>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>1</td>\n",
       "      <td>0.252718</td>\n",
       "      <td>0.049287</td>\n",
       "      <td>0.533373</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "          removed  observed  censored  entrance  at_risk  KM_estimate  \\\n",
       "event_at                                                                \n",
       "0.0             0         0         0        15       15     1.000000   \n",
       "10.0            2         1         1         0       15     0.933333   \n",
       "13.0            1         1         0         0       13     0.861538   \n",
       "18.0            1         1         0         0       12     0.789744   \n",
       "25.0            1         0         1         0       11     0.789744   \n",
       "29.0            1         1         0         0       10     0.710769   \n",
       "30.0            1         1         0         0        9     0.631795   \n",
       "33.0            1         1         0         0        8     0.552821   \n",
       "46.0            1         1         0         0        7     0.473846   \n",
       "50.0            1         0         1         0        6     0.473846   \n",
       "54.0            1         1         0         0        5     0.379077   \n",
       "68.0            1         0         1         0        4     0.379077   \n",
       "71.0            1         1         0         0        3     0.252718   \n",
       "88.0            1         0         1         0        2     0.252718   \n",
       "95.0            1         0         1         0        1     0.252718   \n",
       "\n",
       "          KM_estimate_lower_0.95  KM_estimate_upper_0.95  \n",
       "event_at                                                  \n",
       "0.0                     1.000000                1.000000  \n",
       "10.0                    0.612641                0.990332  \n",
       "13.0                    0.549702                0.963561  \n",
       "18.0                    0.479128                0.927070  \n",
       "25.0                    0.479128                0.927070  \n",
       "29.0                    0.398260                0.881081  \n",
       "30.0                    0.325863                0.828573  \n",
       "33.0                    0.260079                0.770388  \n",
       "46.0                    0.200256                0.706899  \n",
       "50.0                    0.200256                0.706899  \n",
       "54.0                    0.129991                0.630543  \n",
       "68.0                    0.129991                0.630543  \n",
       "71.0                    0.049287                0.533373  \n",
       "88.0                    0.049287                0.533373  \n",
       "95.0                    0.049287                0.533373  "
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# 计算生存率\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "\n",
    "\n",
    "sf2 = kmf2.survival_function_\n",
    "\n",
    "def _kmf_se(kmf: KaplanMeierFitter) -> pd.DataFrame:\n",
    "    \"\"\"Deprecated.\"\"\"\n",
    "    summary = kmf.event_table.merge(kmf.survival_function_, left_index=True, right_index=True).merge(kmf.confidence_interval_, left_index=True, right_index=True)\n",
    "    # 获取生存函数和置信区间\n",
    "    S = kmf.survival_function_\n",
    "    CI_lower = kmf.confidence_interval_.iloc[:, 0]\n",
    "    CI_upper = kmf.confidence_interval_.iloc[:, 1]\n",
    "\n",
    "    # 对数对数变换\n",
    "    loglog_S = -np.log(S)\n",
    "    loglog_CI_lower = np.log(-np.log(CI_lower))\n",
    "    loglog_CI_upper = np.log(-np.log(CI_upper))\n",
    "\n",
    "    # 计算对数对数尺度上的标准误\n",
    "    z_value = 1.96  # 95% 置信水平对应的 Z 分位数\n",
    "    width = loglog_CI_upper - loglog_CI_lower\n",
    "    SE_loglog = width / (2 * z_value)\n",
    "\n",
    "    # 转换回原始尺度的标准误\n",
    "    SE = SE_loglog * (1 / (S * np.log(S)))\n",
    "\n",
    "\n",
    "    temp = pd.DataFrame(SE_loglog).merge( 1 / (S * np.log(S)) , left_index=True, right_index=True)\n",
    "\n",
    "    SE = temp.iloc[:, 0] * temp.iloc[:, 1]\n",
    "\n",
    "    print(temp.iloc[:, 1])\n",
    "\n",
    "    print(SE)\n",
    "\n",
    "    return summary\n",
    "\n",
    "kmf2.event_table.merge(sf2, left_index=True, right_index=True).merge(kmf2.confidence_interval_, left_index=True, right_index=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "乳腺肿瘤直径 $ \\geqslant $ 5cm :"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>removed</th>\n",
       "      <th>observed</th>\n",
       "      <th>censored</th>\n",
       "      <th>entrance</th>\n",
       "      <th>at_risk</th>\n",
       "      <th>KM_estimate</th>\n",
       "      <th>KM_estimate_lower_0.95</th>\n",
       "      <th>KM_estimate_upper_0.95</th>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>event_at</th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0.0</th>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>21</td>\n",
       "      <td>21</td>\n",
       "      <td>1.000000</td>\n",
       "      <td>1.000000</td>\n",
       "      <td>1.000000</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>5.0</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>21</td>\n",
       "      <td>0.952381</td>\n",
       "      <td>0.707207</td>\n",
       "      <td>0.993152</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9.0</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>20</td>\n",
       "      <td>0.904762</td>\n",
       "      <td>0.670046</td>\n",
       "      <td>0.975294</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>13.0</th>\n",
       "      <td>2</td>\n",
       "      <td>2</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>19</td>\n",
       "      <td>0.809524</td>\n",
       "      <td>0.568905</td>\n",
       "      <td>0.923889</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>14.0</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>17</td>\n",
       "      <td>0.761905</td>\n",
       "      <td>0.519391</td>\n",
       "      <td>0.893257</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>15.0</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>16</td>\n",
       "      <td>0.714286</td>\n",
       "      <td>0.471525</td>\n",
       "      <td>0.860197</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>19.0</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>15</td>\n",
       "      <td>0.666667</td>\n",
       "      <td>0.425350</td>\n",
       "      <td>0.825044</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>20.0</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>14</td>\n",
       "      <td>0.619048</td>\n",
       "      <td>0.380844</td>\n",
       "      <td>0.788012</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>21.0</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>13</td>\n",
       "      <td>0.571429</td>\n",
       "      <td>0.337977</td>\n",
       "      <td>0.749241</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>22.0</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>12</td>\n",
       "      <td>0.523810</td>\n",
       "      <td>0.296739</td>\n",
       "      <td>0.708812</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>24.0</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>11</td>\n",
       "      <td>0.476190</td>\n",
       "      <td>0.257141</td>\n",
       "      <td>0.666766</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>25.0</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>10</td>\n",
       "      <td>0.428571</td>\n",
       "      <td>0.219225</td>\n",
       "      <td>0.623104</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>26.0</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>9</td>\n",
       "      <td>0.380952</td>\n",
       "      <td>0.183067</td>\n",
       "      <td>0.577789</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>27.0</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>8</td>\n",
       "      <td>0.333333</td>\n",
       "      <td>0.148785</td>\n",
       "      <td>0.530739</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>28.0</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>7</td>\n",
       "      <td>0.285714</td>\n",
       "      <td>0.116561</td>\n",
       "      <td>0.481820</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>32.0</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>6</td>\n",
       "      <td>0.238095</td>\n",
       "      <td>0.086660</td>\n",
       "      <td>0.430824</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>47.0</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>5</td>\n",
       "      <td>0.190476</td>\n",
       "      <td>0.059482</td>\n",
       "      <td>0.377435</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>52.0</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>4</td>\n",
       "      <td>0.142857</td>\n",
       "      <td>0.035657</td>\n",
       "      <td>0.321162</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>54.0</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>3</td>\n",
       "      <td>0.095238</td>\n",
       "      <td>0.016259</td>\n",
       "      <td>0.261250</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>60.0</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>2</td>\n",
       "      <td>0.047619</td>\n",
       "      <td>0.003324</td>\n",
       "      <td>0.197045</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>86.0</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>1</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "          removed  observed  censored  entrance  at_risk  KM_estimate  \\\n",
       "event_at                                                                \n",
       "0.0             0         0         0        21       21     1.000000   \n",
       "5.0             1         1         0         0       21     0.952381   \n",
       "9.0             1         1         0         0       20     0.904762   \n",
       "13.0            2         2         0         0       19     0.809524   \n",
       "14.0            1         1         0         0       17     0.761905   \n",
       "15.0            1         1         0         0       16     0.714286   \n",
       "19.0            1         1         0         0       15     0.666667   \n",
       "20.0            1         1         0         0       14     0.619048   \n",
       "21.0            1         1         0         0       13     0.571429   \n",
       "22.0            1         1         0         0       12     0.523810   \n",
       "24.0            1         1         0         0       11     0.476190   \n",
       "25.0            1         1         0         0       10     0.428571   \n",
       "26.0            1         1         0         0        9     0.380952   \n",
       "27.0            1         1         0         0        8     0.333333   \n",
       "28.0            1         1         0         0        7     0.285714   \n",
       "32.0            1         1         0         0        6     0.238095   \n",
       "47.0            1         1         0         0        5     0.190476   \n",
       "52.0            1         1         0         0        4     0.142857   \n",
       "54.0            1         1         0         0        3     0.095238   \n",
       "60.0            1         1         0         0        2     0.047619   \n",
       "86.0            1         1         0         0        1     0.000000   \n",
       "\n",
       "          KM_estimate_lower_0.95  KM_estimate_upper_0.95  \n",
       "event_at                                                  \n",
       "0.0                     1.000000                1.000000  \n",
       "5.0                     0.707207                0.993152  \n",
       "9.0                     0.670046                0.975294  \n",
       "13.0                    0.568905                0.923889  \n",
       "14.0                    0.519391                0.893257  \n",
       "15.0                    0.471525                0.860197  \n",
       "19.0                    0.425350                0.825044  \n",
       "20.0                    0.380844                0.788012  \n",
       "21.0                    0.337977                0.749241  \n",
       "22.0                    0.296739                0.708812  \n",
       "24.0                    0.257141                0.666766  \n",
       "25.0                    0.219225                0.623104  \n",
       "26.0                    0.183067                0.577789  \n",
       "27.0                    0.148785                0.530739  \n",
       "28.0                    0.116561                0.481820  \n",
       "32.0                    0.086660                0.430824  \n",
       "47.0                    0.059482                0.377435  \n",
       "52.0                    0.035657                0.321162  \n",
       "54.0                    0.016259                0.261250  \n",
       "60.0                    0.003324                0.197045  \n",
       "86.0                    0.000000                0.000000  "
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "sf5 = kmf5.survival_function_\n",
    "\n",
    "kmf5.event_table.merge(sf5, left_index=True, right_index=True).merge(kmf5.confidence_interval_, left_index=True, right_index=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 生存曲线"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAiMAAAGwCAYAAAB7MGXBAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjkuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8hTgPZAAAACXBIWXMAAA9hAAAPYQGoP6dpAABJTklEQVR4nO3deXhU5d0+8PucMzvZIMuEhGCiQCCySiQGVKBSA+WiIvjKa1EhWlQMKqZUjRXiUo27YItg9WWpyBuUovZl/WlqQBBFobiBUBCMIklMLSRkmZmz/P6IDBmyTjIzZzJzf65rvMhZnzmJmTvPec7zFTRN00BERESkE1HvBhAREVF4YxghIiIiXTGMEBERka4YRoiIiEhXDCNERESkK4YRIiIi0hXDCBEREenKoHcDOkJVVfzwww+IjIyEIAh6N4eIiIg6QNM01NTUICkpCaLYev9HtwgjP/zwA1JSUvRuBhEREXXCd999hz59+rS6vluEkcjISACNbyYqKkrn1hAREVFHVFdXIyUlxf053ppuEUbO3pqJiopiGCEiIupm2htiwQGsREREpCuGESIiItIVwwgRERHpqluMGSEiCjeKosDlcundDKI2GY1GSJLU5eMwjBARBRFN01BeXo5Tp07p3RSiDomJiUFiYmKX5gFjGCEiCiJng0hCQgJsNhsneqSgpWka6urqUFlZCQDo3bt3p4/FMEJEFCQURXEHkdjYWL2bQ9Quq9UKAKisrERCQkKnb9lwACsRUZA4O0bEZrPp3BKijjv789qVMU4MI0REQYa3Zqg78cXPK8MIERER6crrMLJjxw5MmTIFSUlJEAQBb7/9drv7lJaW4pJLLoHZbEa/fv2watWqTjSViIiIQpHXYaS2thbDhg3D0qVLO7T9sWPHMHnyZIwfPx779+/H/Pnz8dvf/hbbtm3zurFERET+MHv2bEydOlXvZoQtr5+mmTRpEiZNmtTh7ZcvX460tDQ899xzAIBBgwZh586deOGFF5CTk+Pt6X1GU1XU19X4/0RGG+DH+79Wo8T7y0Sku9mzZ+PUqVMeveXr16/HjTfeiMcffxxffPEFVq9ejdtvvx3Lly/32DcvLw8vvfQSZs2a5fee8+PHjyMtLQ3//Oc/MXz4cPfyJUuWQNM0v54baPk6dVZRURE2bNiAr7/+GlarFaNHj8ZTTz2F9PT0rjc0wPz+aO/u3bsxYcIEj2U5OTmYP39+q/s4HA44HA7319XV1T5vV31dDWzP9vX5cc/3iToA/+UsBOCfwDA8JQZv3TmagYSIgsqrr76KvLw8LF++HLm5uZg9ezZSUlJQXFyMF154wf1IaENDA9auXYu+ff3/+7gt0dHRup6/M7Zv3468vDxceumlkGUZDz74IK6++mocOHAAPXr00Lt5XvH7ANby8nLY7XaPZXa7HdXV1aivr29xn6KiIkRHR7tfKSkp/m6m31wqHoYVjvY37KT9353CyVMNfjs+EelH0zTUOWVdXl3pJXj66adx1113obi4GLm5ue7ll1xyCVJSUrBhwwb3sg0bNqBv374YMWJEh4+vqiqKioqQlpYGq9WKYcOGYf369e71//nPfzBz5kzEx8fDarWif//+WLlyJQAgLS0NADBixAgIgoBx48YBaH6bZty4cbjrrrswf/589OzZE3a7Ha+88gpqa2uRm5uLyMhI9OvXD1u2bHHvoygKbr31Vne70tPTsWTJEvf6hx9+GKtXr8Y777wDQRAgCAJKS0sBAN999x2uv/56xMTEoFevXrjmmmtw/PjxNq/D1q1bMXv2bFx88cUYNmwYVq1ahbKyMuzdu9e9zalTp3D77bfDbrfDYrFg8ODB2LhxIwBg1apViImJwcaNG5Geng6bzYbrrrsOdXV1WL16NVJTU9GzZ0/cfffdUBSlw9+fzgjKSc8KCgqQn5/v/rq6utrngcRqi0TdgjJU/+ffkB2+/zAX5HokF/8CAPDe7FQoUX26fExNkABj418T9S4FE5d8AABwKWqXj01EwafepSBjkT7j6w48mgObyfuPiPvvvx8vvfQSNm7ciKuuuqrZ+ltuuQUrV67EzJkzAQArVqxAbm6u+0O5I4qKirBmzRosX74c/fv3x44dO3DjjTciPj4eY8eOxcKFC3HgwAFs2bIFcXFxOHLkiPuP3z179mDUqFF47733cPHFF8NkMrV6ntWrV+O+++7Dnj17sG7dOsydOxdvvfUWrr32Wjz44IN44YUXcNNNN6GsrAw2mw2qqqJPnz548803ERsbiw8//BC33XYbevfujeuvvx4LFizAwYMHUV1d7Q5HvXr1gsvlQk5ODrKzs/HBBx/AYDDgj3/8IyZOnIjPP/+8zTY2dfr0afcxgcbQNmnSJNTU1GDNmjW46KKLcODAAY+Jyerq6vDiiy+iuLgYNTU1mDZtGq699lrExMRg8+bN+OabbzB9+nSMGTMGM2bM6PD3yFt+DyOJiYmoqKjwWFZRUYGoqCh3N935zGYzzGazX9sliCJsEdGwRfipa85Z6/5ncv3XgOtY149p6gFcOA4wWlHnlN2L610Kzjjk1vdrgUEUYDF2vbgREdFZW7ZswTvvvIOSkhL84he/aHGbG2+8EQUFBfj2228BALt27UJxcXGHw4jD4cATTzyB9957D9nZ2QCACy+8EDt37sTLL7+MsWPHoqysDCNGjEBmZiYAIDU11b1/fHw8ACA2NhaJiYltnmvYsGF46KGHADT+kfzkk08iLi4Oc+bMAQAsWrQIy5Ytw+eff47LLrsMRqMRjzzyiHv/tLQ07N69G2+88Qauv/56REREwGq1wuFweJx7zZo1UFUVr776qvuW+8qVKxETE4PS0lJcffXV7V4XVVUxf/58jBkzBoMHDwYAvPfee9izZw8OHjyIAQMGuK9VUy6XC8uWLcNFF10EALjuuuvw2muvoaKiAhEREcjIyMD48ePx/vvvd+8wkp2djc2bN3sse/fdd90/RGFBNALmyK4dQ3Y0Bhy1eej44sQpHK7wbjCuzWzA5f3iGEiIgpjVKOHAo/oM9Ld24nfD0KFDUVVVhcLCQowaNQoRERHNtomPj8fkyZOxatUqaJqGyZMnIy4ursPnOHLkCOrq6vDLX/7SY7nT6XTf6pk7dy6mT5+Offv24eqrr8bUqVMxevToTr2fsyRJQmxsLIYMGeJednYIwtnaLACwdOlSrFixAmVlZaivr4fT6fQYKNuSzz77DEeOHEFkpOfnRENDA44ePYoPPvjA48GRl19+2d2zdFZeXh6+/PJL7Ny5071s//796NOnjzuItMRms7mDyNn3lJqa6vG9s9vtHu/RH7wOI2fOnMGRI0fcXx87dgz79+9Hr1690LdvXxQUFODEiRP461//CgC444478Oc//xn33XcfbrnlFvzjH//AG2+8gU2bNvnuXQQ7o7nxqZquUpwtLjZIIiLMxg4fxiErqHPIkFX/jxwnos4TBKFTt0r0kpycjPXr12P8+PGYOHEitmzZ0uwDFmi8VTNv3jwA6PA0EWedOXMGALBp0yYkJyd7rDvboz5p0iR8++232Lx5M959911cddVVyMvLw7PPPuvVuYxGz9+rgiB4LDvbi6GqjbfKi4uLsWDBAjz33HPIzs5GZGQknnnmGXz88cftvqeRI0fi9ddfb7YuPj4eJpMJ+/fvdy87fxzmvHnzsHHjRuzYsQN9+pwbEtDa3Qdv3uPZZWffo794/VP+6aefYvz48e6vz47tOPtI1smTJ1FWVuZen5aWhk2bNuHee+/FkiVL0KdPH7z66qu6PtYbaiySBKvJu79iOM6EiPzhggsuwPbt292BZOvWrc0CycSJE+F0OiEIgtefBRkZGTCbzSgrK8PYsWNb3S4+Ph6zZs3CrFmzcMUVV+D3v/89nn32Wff4C38MyNy1axdGjx6NO++8073s6NGjHtuYTKZm577kkkuwbt06JCQkICoqqsVj9+vXr9kyTdNw11134a233kJpaal7cO5ZQ4cOxffff4/Dhw+32TsSDLwOI+PGjWtzlHVLz4iPGzcO//znP709FRERdUMpKSkoLS3F+PHjkZOTg61bt3qslyQJBw8edP/bG5GRkViwYAHuvfdeqKqKyy+/HKdPn8auXbsQFRWFWbNmYdGiRRg5ciQuvvhiOBwObNy4EYMGDQIAJCQkwGq1YuvWrejTpw8sFovPHuvt378//vrXv2Lbtm1IS0vDa6+9hk8++cQjJKSmpmLbtm04dOgQYmNjER0djZkzZ+KZZ57BNddcg0cffRR9+vTBt99+iw0bNuC+++7z6O1oKi8vD2vXrsU777yDyMhIlJeXA2h8TNlqtWLs2LG48sorMX36dDz//PPo168fvv76awiCgIkTJ/rkPfsKa9N0N656wFEDOM64F0lyLUTXmQ6/JNcZiAofByYi/+nTpw9KS0tRVVWFnJycZvNFRUVFtdoL0J7HHnsMCxcuRFFREQYNGoSJEydi06ZN7g99k8mEgoICDB06FFdeeSUkSUJxcTEAwGAw4MUXX8TLL7+MpKQkXHPNNV17o03cfvvtmDZtGmbMmIGsrCz8+9//9uglAYA5c+YgPT0dmZmZiI+Px65du2Cz2bBjxw707dsX06ZNw6BBg3DrrbeioaGhzWu0bNkynD59GuPGjUPv3r3dr3Xr1rm3+dvf/oZLL70UN9xwAzIyMnDffff5/THdzhC0QEw510XV1dWIjo7G6dOnO/3DG3DOWuCJpMZ/X/9XwNqra8eTHUD1icYnagDUyUDGusb7geuvqECUueO50imrOKOZMPiKaxAR0cWBtUTkMw0NDTh27BjS0tJgsVj0bg5Rh7T1c9vRz+/uMzIq3BnMQFQyoP2caEUNQOO/ZWMEFGPHw4iqNcDQUN/ikzlERESBxjDSnRiazr2iAWh8nFeVrFANXoQRWQXQ8uy3RER6KisrQ0ZGRqvrDxw4oPvU8eR7DCNhrMGlAl5OltZRnFSNiDojKSnJ4zHWltZT6GEYCUOi0DjXyMff/BuK0T8DWTmpGhF1hsFgaPExVgptDCNhyCiJ6GkzQTYboBg7PllaR3FSNSIi8gbDSCDIDsDlgx4Igxn4eca/Lh9KEmExSVD91HPBSdWIiKijGEYCYcMc3xwnIQOY+JTPAgkREVEw4KRn/mK0AX0u9e0xKw809rIQERGFEPaM+IsgADduAA7+H2COAAxdKJQnNwBv3Oi7thEREQURhhF/EoTGcR4GC2DkbIpERMFq9uzZOHXqFN5++229mxKWeJuGiIi6bPbs2Zg6darHsvXr18NiseC5557D7NmzIQgC7rjjjmb75uXlQRAEzJ492+/tPH78OARBaDaXyZIlS1os9OprLV2nznr44YchCILHa+DAgT45dqAxjISABgVokL17BX9FIiLqzl599VXMnDkTy5Ytw+9+9zsAjdV8i4uLUV9/bgbohoYGrF27VvdZVaOjoxETE6NrGzrj4osvxsmTJ92vnTt36t2kTmEYCQG3/EPE9C3o8Ovad63I/zQG3aBGIlF407TGopt6vLrw++Hpp5/GXXfdheLiYuTm5rqXX3LJJUhJScGGDRvcyzZs2IC+fftixIgRHT6+qqooKipCWloarFYrhg0bhvXr17vX/+c//8HMmTMRHx8Pq9WK/v37Y+XKlQDgruw7YsQICIKAcePGAWjeYzFu3DjcddddmD9/Pnr27Am73Y5XXnkFtbW1yM3NRWRkJPr164ctW7a491EUBbfeequ7Xenp6ViyZIl7/cMPP4zVq1fjnXfecfdklJaWAgC+++47XH/99YiJiUGvXr1wzTXX4Pjx4+1eC4PBgMTERPcrLi7OY/2pU6dw++23w263w2KxYPDgwdi4cSMAYNWqVYiJicHGjRuRnp4Om82G6667DnV1dVi9ejVSU1PRs2dP3H333X6v9MsxI92U1QBckiBiX2Xn5vP46rQRDlmDyeTjhhGR77jqzlX/DrQHf3BXCffG/fffj5deegkbN27EVVdd1Wz9LbfcgpUrV2LmzJkAgBUrViA3N9f9odwRRUVFWLNmDZYvX47+/ftjx44duPHGGxEfH4+xY8di4cKFOHDgALZs2YK4uDgcOXLE3RuzZ88ejBo1Cu+99x4uvvhimNr4Jbh69Wrcd9992LNnD9atW4e5c+firbfewrXXXosHH3wQL7zwAm666SaUlZXBZrNBVVX06dMHb775JmJjY/Hhhx/itttuQ+/evXH99ddjwYIFOHjwIKqrq93hqFevXnC5XMjJyUF2djY++OADGAwG/PGPf8TEiRPx+eeft9nGf/3rX0hKSoLFYkF2djaKiorcvUyqqmLSpEmoqanBmjVrcNFFF+HAgQOQpHPzS9XV1eHFF19EcXExampqMG3aNFx77bWIiYnB5s2b8c0332D69OkYM2YMZsyY0eHvkbcYRropQRCwZrIFu4+dgtVogLmDhfIaFGDm//Nz435W7/Rvkg4k1tohat+WLVvwzjvvoKSkBL/4xS9a3ObGG29EQUEBvv32WwDArl27UFxc3OEw4nA48MQTT+C9995DdnY2AODCCy/Ezp078fLLL2Ps2LEoKyvDiBEjkJmZCQBITU117x8fHw8AiI2NRWJiYpvnGjZsGB566CEAQEFBAZ588knExcVhzpzGuaMWLVqEZcuW4fPPP8dll10Go9GIRx55xL1/Wloadu/ejTfeeAPXX389IiIiYLVa4XA4PM69Zs0aqKqKV199FcLP80itXLkSMTExKC0txdVXX91i+7KysrBq1Sqkp6fj5MmTeOSRR3DFFVfgyy+/RGRkJN577z3s2bMHBw8exIABA9zXqimXy4Vly5bhoosuAgBcd911eO2111BRUYGIiAhkZGRg/PjxeP/99xlGqGWCIMAiARYDYO7Ed1JS6iG6moQYwQDV0PWnfiRRQK1Txu6jVV0+VrBgrR3ShdHW2EOh17m9NHToUFRVVaGwsBCjRo1CREREs23i4+MxefJkrFq1CpqmYfLkyc1uLbTlyJEjqKurwy9/+UuP5U6n032rZ+7cuZg+fTr27duHq6++GlOnTsXo0aM79X7OkiQJsbGxGDJkiHuZ3W4HAFRWVrqXLV26FCtWrEBZWRnq6+vhdDoxfPjwNs/z2Wef4ciRI4iMjPRY3tDQgKNHj+KDDz7ApEmT3MtffvllzJw502PZ0KFDkZWVhQsuuABvvPEGbr31Vuzfvx99+vRxB5GW2Gw2dxA5+55SU1M9vnd2u93jPfoDw0gYi6zYA2uTz1bFaEN14pguBxKTQURStA1KiNSmYa0d0o0gdOpWiV6Sk5Oxfv16jB8/HhMnTsSWLVuafcACjbdq5s2bB6Dxw9sbZ86cAQBs2rQJycnJHuvMZjMAYNKkSfj222+xefNmvPvuu7jqqquQl5eHZ5991qtzGc+r3SUIgseys70Yqtp4u7y4uBgLFizAc889h+zsbERGRuKZZ57Bxx9/3O57GjlyJF5//fVm6+Lj42EymTye/jkbgs4XExODAQMG4MiRIwAAq9Xa5fd4dtnZ9+gvDCNhTDVEQPn5J0BQnJBcdYAm++TYpg7eNuouWGuHqGMuuOACbN++3R1Itm7d2iyQTJw4EU6nE4IgICcnx6vjZ2RkwGw2o6ysDGPHjm11u/j4eMyaNQuzZs3CFVdcgd///vd49tln3eMv/DEgc9euXRg9ejTuvPNO97KjR496bGMymZqd+5JLLsG6deuQkJCAqKioFo/dkUrGZ86cwdGjR3HTTTcBaOwt+f7773H48OE2e0eCQWh9YpBXVIPV/dIkjmQlIt9ISUlBaWkpKisrkZOTg+rqao/1kiTh4MGDzQZTdkRkZCQWLFiAe++9F6tXr8bRo0exb98+/OlPf8Lq1asBNI7leOedd3DkyBF89dVX2LhxIwYNGgQASEhIgNVqxdatW1FRUYHTp0/75k0D6N+/Pz799FNs27YNhw8fxsKFC/HJJ594bJOamorPP/8chw4dQlVVFVwuF2bOnIm4uDhcc801+OCDD3Ds2DGUlpbi7rvvxvfff9/q+RYsWIDt27fj+PHj+PDDD3HttddCkiTccMMNAICxY8fiyiuvxPTp0/Huu+/i2LFj2LJlC7Zu3eqz9+wrDCNERORzffr0QWlpKaqqqloMJFFRUa32ArTnsccew8KFC1FUVIRBgwZh4sSJ2LRpk/uxXZPJhIKCAgwdOhRXXnklJElCcXExgMZHYV988UW8/PLLSEpKwjXXXNO1N9rE7bffjmnTpmHGjBnIysrCv//9b49eEgCYM2cO0tPTkZmZifj4eOzatQs2mw07duxA3759MW3aNAwaNAi33norGhoa2rxG33//PW644Qakp6fj+uuvR2xsLD766CP3IF0A+Nvf/oZLL70UN9xwAzIyMnDffff5/THdzhC0bjDZRHV1NaKjo3H69OlO//DqwlEDHN4GmCM7NRjMzdUArL2u8d+/We+eWr7epWDft/+B1eTF0zRy41wjAPC3SY2DXwFAlOshuc7gP33GQzU2H3QWzuqdCs44XPjFIDsiOjNSmKiDGhoacOzYMaSlpcFiYQkJ6h7a+rnt6Oc3f7N2N3JDk38rEBUHREWGIHQsjIgKYAVQDzMAwS9NJCIi8gbDSHfTpHqvFUB2Jw5x0AJ8og4AtEIwkBBRMCkrK0NGRkar6w8cOKD71PHkewwj3YHBDCRkAJUHfHbIS8XD+KfqAMCuYCIKHklJSc2K2J2/nkIPw0h3IAjAxKcA2eGxuF5WsL/sFKwmqcOP0rocDgzbNdcfrSQi6jKDwdChx1gptDCMdBeC4B64eo4CVTJDlQzQpI6FEZUTiBIFPX9PMEXkS774eWUYIQ+SfK60t6+mhyeijjGZTBBFET/88IN75s2zs3wSBRtN0+B0OvHjjz9CFMU2C/q1h2GEAACaIEGUaxFV/pF7ma+mhyeijhFFEWlpaTh58iR++EGnmjREXrLZbOjbty9EsfNTlzGMEABAk0xw2pIgaI2T4fh6engi6hiTyYS+fftCluWgnJyKqClJkmAwGLrcg8cwQm6aZMLZGfBEAFCdOraGKHydLVZ2fsEyolDFMELUQfVOz79SDaIAi5EjgomIuophhKgdkiig1ilj99Eqj+U2swGX94tjICEi6iKGEaJ2mAwikqJtUNRzZZwcsoI6hwxZDfrSTkREQY9hhKgDWppUzqVwLggiIl/o/HM4RERERD7AnpEw1qAAaitP7opyY4XfBpcKFd4/Xmg2iJysiYiIOoRhJIzdUgLUt7rW+vOrc8X50u2RKJySwUBCRETtYhgJM6YAPfhxqKIGrvoaWIw/3wnk1PJERNQKhpEw07SjYu3VrRfOExQnjPUnoRpsXh2/QREwY0ccACDmh+2w/nx8Ti1PREStYRgJY2YDoLXWU2IwQTD0hkHzbrxI0zEoqiECioFTyxMRUdsYRqhVTaeH76imD7uqBitUA6eWJyKitvHRXiIiItIVwwgRERHpimGEiIiIdMUwQkRERLpiGCEiIiJdMYwQERGRrhhGiIiISFecZ4SoC+qd3hcRNIgCLMYAzctPRNQNMIxQwEjyeWX5unG9GkkUUOuUsftoldf72swGXN4vjoGEiOhnDCPkd5ogQZRrEVX+kcfy7lyvxmQQkRRtg6J6N0etQ1ZQ55Ahe7kfEVEoYxghv9MkE5y2JAhN6tyEQr0ak6FzQ65citr+RkREYYRhhALi/Do3rFdDRERnMYyEMVFxeBS200QzIAi6tYeIiMITw0gYG7h9rsfXtTEDcDyzkIGEiIgCqlM3vZcuXYrU1FRYLBZkZWVhz549bW6/ePFipKenw2q1IiUlBffeey8aGho61WDqGk00ozZmQIvrepw6DEF1BLhFREQU7rzuGVm3bh3y8/OxfPlyZGVlYfHixcjJycGhQ4eQkJDQbPu1a9figQcewIoVKzB69GgcPnwYs2fPhiAIeP75533yJsgLgoDjmYUeoUNUHM16SYiIiALF656R559/HnPmzEFubi4yMjKwfPly2Gw2rFixosXtP/zwQ4wZMwa/+c1vkJqaiquvvho33HBDm70pDocD1dXVHi/yIUGAJlncL1Uy690iIiIKY16FEafTib1792LChAnnDiCKmDBhAnbv3t3iPqNHj8bevXvd4eObb77B5s2b8atf/arV8xQVFSE6Otr9SklJ8aaZRERE1I14dZumqqoKiqLAbrd7LLfb7fj6669b3Oc3v/kNqqqqcPnll0PTNMiyjDvuuAMPPvhgq+cpKChAfn6+++vq6moGEiIiohDl96dpSktL8cQTT+Cll15CVlYWjhw5gnvuuQePPfYYFi5c2OI+ZrMZZjNvHXR3DW2UbRFlQFSABpcKFd7Xdwk0s0GEwKeMiIj8wqswEhcXB0mSUFFR4bG8oqICiYmJLe6zcOFC3HTTTfjtb38LABgyZAhqa2tx22234Q9/+ANEkYWDQ9XM/9fWWuvPrwOBaUwXpdsjUTglg4GEiMgPvEoCJpMJI0eORElJiXuZqqooKSlBdnZ2i/vU1dU1CxyS1FggTNNYnyPUmCUgo6ferfC9QxU1cMicxp2IyB+8vk2Tn5+PWbNmITMzE6NGjcLixYtRW1uL3NxcAMDNN9+M5ORkFBUVAQCmTJmC559/HiNGjHDfplm4cCGmTJniDiUUOgQBeHoM4OjAnRdBcTarVyNoTpxKGgvV2MOPrew4h6zijjV79W4GEVFI8zqMzJgxAz/++CMWLVqE8vJyDB8+HFu3bnUPai0rK/PoCXnooYcgCAIeeughnDhxAvHx8ZgyZQoef/xx370LCiqCAFg68pNlMHl8KcqA5HKiwShCNTKoEhGFi04NYJ03bx7mzZvX4rrS0lLPExgMKCwsRGFhYWdORURERCGOo0eJiIhIVyyUR6SDemfwP87sawZRgIW334ioBQwjRAEkiQJqnTJ2H63SuykBZzMbcHm/OAYSImqGYYQ8iIoDZx9g1URz42hU8hmTQURStA2KGl6PtTtkBXUOGXKYvW8i6hiGEfLQtHpvbcwAHM8sZCDxMZMhPIdquRTO00JELQvP34rkQRPNqI0Z0Gx5j1OHIagOHVpEREThhD0jBAgCjmcWuoOHqDg8ekiIiIj8iWGEGgkCNMkCAGBnOhERBRLDCAUdSa73XCAYoBos+jSGiIj8jmGEgoYmSBDlWkSVf+SxXDHaUJ04hoGEiChEMYxQ0NAkE5y2pGbF8yRXHaDJOraMiIj8iWGEgoommdB0JgoRAFSnTq0hIqJA4KO9REREpCuGESIiItIVwwgRERHpimGEiIiIdMUwQkRERLpiGCEiIiJdMYwQERGRrjjPCLVJVBw+qVWjiWZAEHxwJCIiCjUMI9QmX1XvrY0ZgOOZhd06kDhkz1hmNogQuvH7ISIKFgwj1IwmmlEbMwA9Th322TF7nDoMQXW4KwN3R3es2evxdbo9EoVTMhhIiIi6iGEkEGSHb44jSIDB7JtjtXkeAcczCyGoXW+3qDh80rvSrJLv+fxU2ddsEJFuj8Shippm6w5V1MAhq7AYJZ+fN1TVO5X2Nwowgyjwe0ikM4YRfxINgKkH4KwFFB/UV3HWAlHJAQskvujF6Op4k9Yq+Z7PX5V9BUFA4ZQMj1s0Dllt1ktCbZNEAbVOGbuPVundlGZsZgMu7xfHQEKkI4YRfzJagQvHAaoPKs666oFjOwAt+P6y9KeWKvmez9+VfQWBfzl3lckgIinaBkXV2t84gByygjqHDDnI2kUUbhhG/M1o1bsF3d75lXzPx8q+3YPJEJwzCbgUXzwvRkRdEZy/HYiIiChsMIwQERGRrhhGiIiISFcMI0RERKQrhhEiIiLSFcMIERER6YphhIiIiHTFMEJERES6YhghIiIiXXEG1u6madE9lwJJrocoGiB6mSs1QYImmXzcuLaJisOjVo0mmgFWvCUiCnsMI91FS0X3XCoMci0MggRJ8y6MiHItnLakgAaS86v31sYMwPHMQp8FkmaVff1UyZeIiHyLYaS7aKHonuKQUaH8iAizARZTxwu5SXI9oso/gqApbdZ88QVNNKM2ZgB6nDrcbF2PU4chqI4uVwdurbKvvyr5EhGRbzGMdCfNiu7JUIx1UIxGqMFaVVYQcDyzEIJ67vaSqDia9ZJ0RUuVff1dyZeIiHyHYYT8TxA8ej/8USP1/Mq+rORLRNR98GkaIiIi0hXDCBEREemKYYSIiIh0xTBCREREumIYISIiIl0xjBAREZGuGEaIiIhIV5xnhKgLHLL3s6aYDSIE1uQhInJjGCHqgjvW7PV6n3R7JAqnZDCQBJF6p9L+RtQigyjAEqwzQFO3wTBC5CWzQUS6PRKHKmo6tf+hiho4ZJW/wIOAJAqodcrYfbRK76Z0WzazAZf3i+PPM3UJwwiRlwRBQOGUDK9v0ThktVM9KeQ/JoOIpGgbFNXfJSNDk0NWUOeQIfP6URcxjJCuRMXhUatGE82AD29fSHK95wLB4JMqvoLArulQYTJwHH9XuBR/VJuicMMwQro6v3pvbcwAHM8s7HIg0QQJolyLqPKPPJYrRhuqE8f4JJAQEZFv8E8CCjhNNKM2ZkCL63qcOgxBdXT9HJIJTlsSFGOE+6WKJkiuOkCTu3x8IiLyHfaMhACH7N2TAJJLgayoEIx+alB7BAHHMws9QoeoOJr1knSVJpnQ9E62CACq06fnICKirmMY6cYMogCb2YA6h+zVfVvJJcNQ50SkUYWk10+AIECTzt0q4V1nIqLw1amPoqVLl+KZZ55BeXk5hg0bhj/96U8YNWpUq9ufOnUKf/jDH7Bhwwb89NNPuOCCC7B48WL86le/6nTDCbAYJVzeL87rkewNZ07jeLkEVQM4BJOIiPTmdRhZt24d8vPzsXz5cmRlZWHx4sXIycnBoUOHkJCQ0Gx7p9OJX/7yl0hISMD69euRnJyMb7/9FjExMb5of9jr1BMdDg4VIiKi4OF1GHn++ecxZ84c5ObmAgCWL1+OTZs2YcWKFXjggQeabb9ixQr89NNP+PDDD2E0Ng5SSE1NbfMcDocDDse58QTV1dXeNpOIiIi6Ca/+RHY6ndi7dy8mTJhw7gCiiAkTJmD37t0t7vP3v/8d2dnZyMvLg91ux+DBg/HEE09AUVofdFlUVITo6Gj3KyUlxZtmEhERUTfiVRipqqqCoiiw2+0ey+12O8rLy1vc55tvvsH69euhKAo2b96MhQsX4rnnnsMf//jHVs9TUFCA06dPu1/fffedN80kIiKibsTvz1KoqoqEhAT85S9/gSRJGDlyJE6cOIFnnnkGhYWFLe5jNpthNpv93TQiIiIKAl6Fkbi4OEiShIqKCo/lFRUVSExMbHGf3r17w2g0QpLODbQcNGgQysvL4XQ6YTKZOtFsIiIiChVe3aYxmUwYOXIkSkpK3MtUVUVJSQmys7Nb3GfMmDE4cuQIVPXcTBKHDx9G7969GUSIiIjI++ng8/Pz8corr2D16tU4ePAg5s6di9raWvfTNTfffDMKCgrc28+dOxc//fQT7rnnHhw+fBibNm3CE088gby8PN+9CyIiIuq2vB4zMmPGDPz4449YtGgRysvLMXz4cGzdutU9qLWsrAyieC7jpKSkYNu2bbj33nsxdOhQJCcn45577sH999/vu3dBRERE3VanBrDOmzcP8+bNa3FdaWlps2XZ2dn46KOPmm9MpANJrj/3hWBgBV8iIp2xNg0FHVFxeNSq0UQzIAhdPq4mSBDlWkSVnwvGitGG6sQxDCRERDpiGKGgc3713tqYATieWdjlQKJJJjhtSRC0xgn3BMUJyVUHaHKXjktERF3DIiVhTFSdEOV690tQnLq1RRPNqI0Z0OK6HqcOQ1AdLa7z+jySCarBCtVghSbxaS4iomDAnpFwJBogG6ywKE5IrnO9AqJcC6ctSZ8PaUHA8cxCj9AhKo5mvSRERBR6GEbCkdGKqoRsOIwCLKbGyegkuR5R5R9B0BRoerVLEKBJ58ZuqG1sSkREoYNhJEypkgWK0QjVKLW/MRERkR8xjBDpwCGHX7+P2SBC8MFTUUQUehhGiHRwx5q9ejch4NLtkSicksFAQkTNMIyEMYesuP8tuRQ4ZRWyoEL1wWgNUQCMEh/WaspsEJFuj8Shihq9m6KLQxU1cMgqLLw1SETnYRgJQwZRgM1sQJ1DhktpDB6SS0YPlwJZk6GoXZ93o0FWENvDzEDShCAIKJySEXa3aByyGpY9QUTUcQwjYchilHB5vzjIapPnZhw1kKRowBwBGG1dOr5DVvHVidNQdXssJ3gJgsCeASKi8zCMhKnmH4gGwCgCRqnxRUREFCDsQyciIiJdMYwQERGRrnibhoiIuqTeqbS/EQU1g6jveDaGEQp7klzvuUAwQDVYWt6YiNwkUUCtU8buo1V6N4W6yGY24PJ+cboFEoYR6hZExdHm7CeaaAa8nExLEySIci2iyj/yWK4YbahOHMNAQtQOk0FEUrQNCh+d69YcsoI6h+z5hGWAMYxQt9Be9d7amAE4nlnoVSDRJBOctiQI2rkuZkFxQnLVAVrX51ohCgcmA4cehoKzc07phT9F5El2AK66xpfs0LUpmmhGbcyADm3b49RhCKr37dUkE1SD1f3SJJPXxyAioq5hzwg1Eg2AqQfgrAUUZ+MyZy0QlQwYzPq0SRBwPLOwzZAhKo52e02IiCi4MYxQI6MVuHAccHYqeFc9cGwHoOk8Sl4QoEmtj90Ir4nViYhCE8MInWO06t0CIiIKQxwzQkRERLpiGCEiIiJdMYwQERGRrhhGiIiISFcMI0RERKQrhhEiIiLSFcMIERER6YphhIiIiHTFMEJERES64gysRC2Q5HrfHEgwQDW0Pp09ERExjJAfuRT/17URm5S9lh31UKVzX6uiGRAEr44nKAIsDTWwfP+h53kEwCB535GoGG2oThzDQEJE1AaGEfI5SRRgMUlocCqQFdmv5xKbBJ6hu/I81lVH9scXQx70MpCIaDAkQDivQGCDrKCX0QyjF4FEUJyQXHWA5t9r0J045OArbWg2iBC8DK1E5FsMI9Q22eH1LiYAQxKMUNQmP16CBBjMvmvXWZoG5eggSD8ebLYqquZfuKSPDehir4RDVvHVidOQJQMkQ8fDiAgAqrNL5w41d6zZq3cTmkm3R6JwSgYDCZGOGEaoZaIBMPUAnLWA4v0Hqun8Bc5aICrZP4Fk0tOeoUluAN64EQBgNUiAUfL9OanDzAYR6fZIHKqo0bspLTpUUQOHrMLCnxMi3TCMUMuMVuDCcYDqg1sMrnrg2A5A89MYEkEAjByTEawEQUDhlIygu0XjkNWg7KkhCkcMI9Q6o1XvFlCIEASBPQ9E1CrOM0JERES6YhghIiIiXTGMEBERka4YRoiIiEhXDCNERESkK4YRIiIi0hXDCBEREemKYYSIiIh0xTBCREREumIYISIiIl1xOngKbXKD59cGc2MtmwCS5PqAnq/bEAxQu1hRmYhCA8MIhbafq/e6JWQAE58KSCDRBAmiXIuo8o/8fq7uSDHaUJ04hoGEiBhGKIBkh+fXgtTYU+FrBnNj6Kg80Hxd5YHGdgSgyq8mmeC0JUHwV7XibkxQnJBcdYDmg6rQRNTtMYyQ/4kGwNQDcNYCivPccmctEJXs+0AiCI29H03Dj9zQvJckADTJBC3gZw1+IgCozvY2I6IwwTBC/me0AheOA9QmfwW76oFjOwB/9RoIQkB6P4iIqOsYRigwjFa9W0BEREGKj/YSERGRrhhGiIiISFedCiNLly5FamoqLBYLsrKysGfPng7tV1xcDEEQMHXq1M6cloiIiEKQ12Fk3bp1yM/PR2FhIfbt24dhw4YhJycHlZWVbe53/PhxLFiwAFdccUWnG0tEREShx+sw8vzzz2POnDnIzc1FRkYGli9fDpvNhhUrVrS6j6IomDlzJh555BFceOGFXWowERERhRavwojT6cTevXsxYcKEcwcQRUyYMAG7d+9udb9HH30UCQkJuPXWWzt0HofDgerqao8Xkd5cigKHrLpfLkXVu0lERCHBqzBSVVUFRVFgt9s9ltvtdpSXl7e4z86dO/E///M/eOWVVzp8nqKiIkRHR7tfKSkp3jSTyKckUYDFJEFWNNQ7Zffr37UOBhIiIh/w6zwjNTU1uOmmm/DKK68gLi6uw/sVFBQgPz/f/XV1dTUDCenGJIkYkhwNRT03l6pDVvHVidNQOb0qEVGXeRVG4uLiIEkSKioqPJZXVFQgMTGx2fZHjx7F8ePHMWXKFPcyVW38S9JgMODQoUO46KKLmu1nNpthNvuhZglRJ5kkEZD0bgURUWjy6jaNyWTCyJEjUVJS4l6mqipKSkqQnZ3dbPuBAwfiiy++wP79+92vX//61xg/fjz279/P3g4iIiLy/jZNfn4+Zs2ahczMTIwaNQqLFy9GbW0tcnNzAQA333wzkpOTUVRUBIvFgsGDB3vsHxMTAwDNlhMFjNxw7t8Gc2MdGyIi0o3XYWTGjBn48ccfsWjRIpSXl2P48OHYunWre1BrWVkZRJETu1IHNa2sCwCC5PsqvudrWr03IaOxwi8DSVhzyByI3FlmgwiB//9QF3VqAOu8efMwb968FteVlpa2ue+qVas6c0oKNaIBMPUAnLWA0qSUvLMWiEr2fSAxmBuDR+UBz+WVBxoDESv86kKS6/13cMEA1dCx7+sda/b6rx0hLt0eicIpGQwk1CWs2kv6MFqBC8cBqnxumaseOLYD0BTfn08QGntAzvbEyA2ePSQUUJogQZRrEVX+kd/OoRhtqE4c02ogMRtEpNsjcaiixm9tCAeHKmrgkFVYjBzhTZ3HMEL6MVoDez5BYA9IkNAkE5y2JAj+CJ4ABMUJyVUHaHLr2wgCCqdk8BZNJzlklT1K5DMMI0SkC00ywV/TtIgAoDrb2wyCIPAveqIgwJGmREREpCuGESIiItIVwwgRERHpimGEiIiIdMUwQkRERLpiGCEiIiJdMYwQERGRrhhGiIiISFcMI0RERKQrzsBKwUePSr5ERKQbhhEKHoGu5HuW3OD5tcHcWMeGiIgCgmGEgkegK/medX713oSMxgq/DCRERAHBMSMUXIxWwBx57uWvyr4Gc2PoaEnlgea3ioiIyG/YM0LhSRAaez+ahg65oXkvSTtcivc9NqIAGCX+HUBEdBbDCIUvQQCMlk7tKokCLCYJDU4FsiK3v0MTDbKC2B5mBhIiop8xjBB1gkkSMSQ5GoqqebWfQ1bx1YnT8HI36iRJrm9/I8EA1dC5UEpEvsEwQtRJJkkEJL1bQS3RBAmiXIuo8o/a3VYx2lCdOIaBhEhHDCNEFHI0yQSnLQlCO09hCYoTkqsO0Ly71UZEvsUwQkQhSZNMaO9umAgAqrOdrYjI3ziCjoiIiHTFMEJERES6YhghIiIiXTGMEBERka44gJW6h/amZ2dlXyKibothhIJba5V8z+fvyr5EROQ3DCMU3Fqq5Hu+QFT2JSIiv2EYoeDnr8q9bZEbfHMcg7mxBg4REbWKYYSoJV5W721VQkZjdWAGEiKiVvFpGqKzDObG8OBLlQfaH3xLRBTm2DNCdJYgNPZi+CI8yA2+610hIgpxDCNETQkCYGT1ViKiQOJtGiIiItIVwwgRERHpimGEiIiIdMUwQkRERLriAFYiCnuSXB+4kwkGqAYOkiZqimGESAcuJfymrhcFwCgFV2esJkgQ5VpElX8UsHMqRhuqE8cwkBA1wTBCoeP8+UGCsJKvJAqwmCQ0OBXIShv1dkJQg6wgtoc5qAKJJpngtCVBCFBdI0FxQnLVAVpofe8dsqp3E6gLHC4FDlmBpmm6tYFhhLq/1ir7BmElX5MkYkhyNBRVv//p9eCQVXx14jSC8W1rkgmBapYIAGob1ae7qTvW7NW7CeQDY9MTEGkx6nJuhhHq/lqq7BvElXxNkghIereCqGvMBhHp9kgcqqjRuykUAhhGKDToUdmXKIwJgoDCKRm8RRMCGpwKzjhdsBr1+yuJYYSIiDpFEARYdPwAI9/QNMClqhB0rC7OMELkb3KD3i3Qn6xAVBwQFRmC0DiAVRPNjbWAiCjsMYwQ+Rur98IKIPu8ZbUxA3A8s5CBhIg4AyuRXxjMQEKG3q0Iaj1OHYagOtrfkIhCHntGiPxBEICJTzWf+yRM1csK9pedgtUkwSK4MHD7XL2bRERBhGGEyF8EATByls1GClTJDFUyQBXYIUtEnvhbgYiIiHTFMEJERES64m0aCm3doF4NEVG4Yxih0NSN6tUQEYU7hhEKTd2sXg0RUThjGKHQxXo1RETdQqcGsC5duhSpqamwWCzIysrCnj17Wt32lVdewRVXXIGePXuiZ8+emDBhQpvbExERUXjxOoysW7cO+fn5KCwsxL59+zBs2DDk5OSgsrKyxe1LS0txww034P3338fu3buRkpKCq6++GidOnOhy44mIiKj78zqMPP/885gzZw5yc3ORkZGB5cuXw2azYcWKFS1u//rrr+POO+/E8OHDMXDgQLz66qtQVRUlJSVdbjwRERF1f16NGXE6ndi7dy8KCgrcy0RRxIQJE7B79+4OHaOurg4ulwu9evVqdRuHwwGH49wjmdXV1d40k4iClEtRIEJ1f+2UVaia2sYe/icKgFHilEtEevIqjFRVVUFRFNjtdo/ldrsdX3/9dYeOcf/99yMpKQkTJkxodZuioiI88sgj3jSNiIKYJAqwmCQ0OBWoyrmnmeqdClRJbmNP/2uQFcT2MAc0kEhyfcDORToSDFANLAnREQF9mubJJ59EcXExSktLYbG0/g0qKChAfn6+++vq6mqkpKQEoolE5AcmScSQ5GgoqgbIDcBHjcuHJ1kAf/2ylsyN9YHa4JBVfHXiNFTNP004nyZIEOVaRJV/FJgTkq4Uow3ViWMYSDrAqzASFxcHSZJQUVHhsbyiogKJiYlt7vvss8/iySefxHvvvYehQ4e2ua3ZbIbZzEmpiEKJSRIBCfj5PwAA64ZZ/jthQkZj5eR2AkkgaZIJTlsSBM51E/IExQnJVQdo+vb8dRdehRGTyYSRI0eipKQEU6dOBQD3YNR58+a1ut/TTz+Nxx9/HNu2bUNmZmaXGkxE3ZzB3BgUKg/49zyVBxrLAQRZ5WRNMiFAHTGkIxEAVGd7m9HPvL5Nk5+fj1mzZiEzMxOjRo3C4sWLUVtbi9zcXADAzTffjOTkZBQVFQEAnnrqKSxatAhr165FamoqysvLAQARERGIiIjw4Vsh6qCm9WpYqybwBKGxx+L8ukG+IjcAb9zon2MTkV94HUZmzJiBH3/8EYsWLUJ5eTmGDx+OrVu3uge1lpWVQRTPDQRbtmwZnE4nrrvuOo/jFBYW4uGHH+5a64m80VK9Gtaq0YcgBF2PBRHpp1MDWOfNm9fqbZnS0lKPr48fP96ZUxD53vn1alirhogoKLA2DYUX1qshIgo6nOmHiIiIdMUwQkRERLpiGCEiIiJdMYwQERGRrhhGiIiISFcMI0RERKQrhhEiIiLSFcMIERER6YqTnhGdXyOF9WqIiAKKYYTCV0u1agDWqwklckM76xVAYw1dIr0xjFD4Or9WDcB6NaGmneq9VgBDIvvjyMhFgWkPEbWIY0YovBmtgDny3Iu1a7o/gxlIyOjw5lE1/4KoOtrfkIj8hj0jRBRaBAGY+FTzsUDnkxva7TkhosBgGCGi0CMIgNGidyuIqIMYRogo7LkUFaqs6t2MbkkUAKPEO/7UNQwjRBT2ZEWD0ym3vyE10yAriO1hZiChLmEYIaKwN7xvDGDgbR1vOWQVX504DZVPR1MXMYwQUdizGiTAKOndDKKwxX41IiIi0hV7Roha0t5joR3FqeWJiNrFMELUVGtTxHcWp5YnImoXwwhRUy1NEd9ZnFqeiKhDGEaIzscp4YmIAooDWImIiEhXDCNERESkK96mISKSGwJ3LoO5sXYOEbkxjBARBbJ6b0JGY1VhBhIiN96mIaLwZDA3BoNAqzzgu3lsiEIEe0aIKDwJQmMPRaCCgdwQ2B4Yom6EYYTI3/hXcMuCYXZaQQCMLJBHpDeGESJ/8fVsrqGGs9MS0c8YRoj8xZezuYYazk5LRE0wjBD5E2dzJSJqF5+mISIiIl0xjBAREZGuGEaIiIhIVwwjREREpCuGESIiItIVn6YhIv34c0K4YJhUjYg6hGGEiAIvEBPCcVI1om6DYYSIAs/fE8JxUjWiboVhhIj0Ec4TwskNerfAN2QFouKA4nTCpXAIYlOi3ABFUdDgVKAEeSh2yPq3j2GEiCjQQqR6rxVAtt6NCGI1PVLxbezlUFSX3k1pl81sgEEUdDs/wwgRUSAYzEBCBlB5QO+WUIBE1h7H2IsigAi73k1pl0EUYDFK+p1ftzMTEYUTQQAmPuXfJ4goOMgN7t6vCJMBMPOjtj28QkREgSIIgNGidyuIgg7DCBGFro70QnA+EiLdMYwQUejxZh4TzkdCpDuGESIKPR2dx4TzkRAFBYYRIgpN4TyPCVE3w1lqiIiISFcMI0RERKQr3qYhIgrk3B98eoeoGYYRIgpfgagefD4+vUPUDMMIEYUvf1cPPh+f3iFqUafGjCxduhSpqamwWCzIysrCnj172tz+zTffxMCBA2GxWDBkyBBs3ry5U40lIvI5oxUwRwbmxSd8iFrkdRhZt24d8vPzUVhYiH379mHYsGHIyclBZWVli9t/+OGHuOGGG3Drrbfin//8J6ZOnYqpU6fiyy+/7HLjiYiIqPsTNE3TvNkhKysLl156Kf785z8DAFRVRUpKCu666y488MADzbafMWMGamtrsXHjRveyyy67DMOHD8fy5cs7dM7q6mpER0fj9OnTiIqK8qa5RETBw1EDHN4GSCaOGQllcgPwxs2N//7dYSAy+Kv2+ktHP7+9GjPidDqxd+9eFBQUuJeJoogJEyZg9+7dLe6ze/du5OfneyzLycnB22+/3ep5HA4HHI5zo9urq6u9aSYRUXDSY8AsBV7Tp7NcDfq1oxvxKoxUVVVBURTY7Z4pz2634+uvv25xn/Ly8ha3Ly8vb/U8RUVFeOSRR7xpGhFR8Av0gFnSh7Pu3L8lPifSEUF5lQoKCjx6U6qrq5GSkqJji4iIfISDWEOfKQJ48AfAWQ/0iNW7Nd2CV2EkLi4OkiShoqLCY3lFRQUSExNb3CcxMdGr7QHAbDbDbOb9VCIi6oYEofF2nKmH3i3pNrx6msZkMmHkyJEoKSlxL1NVFSUlJcjOzm5xn+zsbI/tAeDdd99tdXsiIiIKL17fpsnPz8esWbOQmZmJUaNGYfHixaitrUVubi4A4Oabb0ZycjKKiooAAPfccw/Gjh2L5557DpMnT0ZxcTE+/fRT/OUvf/HtOyEiIqJuyeswMmPGDPz4449YtGgRysvLMXz4cGzdutU9SLWsrAyieK7DZfTo0Vi7di0eeughPPjgg+jfvz/efvttDB482HfvgoiIiLotr+cZ0QPnGSEiIup+Ovr53anp4ImIiIh8hWGEiIiIdMUwQkRERLpiGCEiIiJdMYwQERGRrhhGiIiISFcMI0RERKQrhhEiIiLSFcMIERER6crr6eD1cHaS2Orqap1bQkRERB119nO7vcneu0UYqampAQCkpKTo3BIiIiLyVk1NDaKjo1td3y1q06iqih9++AGRkZEQBMFnx62urkZKSgq+++471rwJMF57/fDa64fXXj+89vrQNA01NTVISkryKKJ7vm7RMyKKIvr06eO340dFRfGHUye89vrhtdcPr71+eO0Dr60ekbM4gJWIiIh0xTBCREREugrrMGI2m1FYWAiz2ax3U8IOr71+eO31w2uvH1774NYtBrASERFR6ArrnhEiIiLSH8MIERER6YphhIiIiHTFMEJERES6CuswsnTpUqSmpsJisSArKwt79uzRu0khpaioCJdeeikiIyORkJCAqVOn4tChQx7bNDQ0IC8vD7GxsYiIiMD06dNRUVGhU4tD15NPPglBEDB//nz3Ml57/zlx4gRuvPFGxMbGwmq1YsiQIfj000/d6zVNw6JFi9C7d29YrVZMmDAB//rXv3RscWhQFAULFy5EWloarFYrLrroIjz22GMedVF47YOUFqaKi4s1k8mkrVixQvvqq6+0OXPmaDExMVpFRYXeTQsZOTk52sqVK7Uvv/xS279/v/arX/1K69u3r3bmzBn3NnfccYeWkpKilZSUaJ9++ql22WWXaaNHj9ax1aFnz549WmpqqjZ06FDtnnvucS/ntfePn376Sbvgggu02bNnax9//LH2zTffaNu2bdOOHDni3ubJJ5/UoqOjtbffflv77LPPtF//+tdaWlqaVl9fr2PLu7/HH39ci42N1TZu3KgdO3ZMe/PNN7WIiAhtyZIl7m147YNT2IaRUaNGaXl5ee6vFUXRkpKStKKiIh1bFdoqKys1ANr27ds1TdO0U6dOaUajUXvzzTfd2xw8eFADoO3evVuvZoaUmpoarX///tq7776rjR071h1GeO395/7779cuv/zyVterqqolJiZqzzzzjHvZqVOnNLPZrP3v//5vIJoYsiZPnqzdcsstHsumTZumzZw5U9M0XvtgFpa3aZxOJ/bu3YsJEya4l4miiAkTJmD37t06tiy0nT59GgDQq1cvAMDevXvhcrk8vg8DBw5E3759+X3wkby8PEyePNnjGgO89v7097//HZmZmfiv//ovJCQkYMSIEXjllVfc648dO4by8nKPax8dHY2srCxe+y4aPXo0SkpKcPjwYQDAZ599hp07d2LSpEkAeO2DWbcolOdrVVVVUBQFdrvdY7ndbsfXX3+tU6tCm6qqmD9/PsaMGYPBgwcDAMrLy2EymRATE+Oxrd1uR3l5uQ6tDC3FxcXYt28fPvnkk2breO3955tvvsGyZcuQn5+PBx98EJ988gnuvvtumEwmzJo1y319W/r9w2vfNQ888ACqq6sxcOBASJIERVHw+OOPY+bMmQDAax/EwjKMUODl5eXhyy+/xM6dO/VuSlj47rvvcM899+Ddd9+FxWLRuzlhRVVVZGZm4oknngAAjBgxAl9++SWWL1+OWbNm6dy60PbGG2/g9ddfx9q1a3HxxRdj//79mD9/PpKSknjtg1xY3qaJi4uDJEnNnhyoqKhAYmKiTq0KXfPmzcPGjRvx/vvvo0+fPu7liYmJcDqdOHXqlMf2/D503d69e1FZWYlLLrkEBoMBBoMB27dvx4svvgiDwQC73c5r7ye9e/dGRkaGx7JBgwahrKwMANzXl79/fO/3v/89HnjgAfz3f/83hgwZgptuugn33nsvioqKAPDaB7OwDCMmkwkjR45ESUmJe5mqqigpKUF2draOLQstmqZh3rx5eOutt/CPf/wDaWlpHutHjhwJo9Ho8X04dOgQysrK+H3ooquuugpffPEF9u/f735lZmZi5syZ7n/z2vvHmDFjmj3CfvjwYVxwwQUAgLS0NCQmJnpc++rqanz88ce89l1UV1cHUfT8WJMkCaqqAuC1D2p6j6DVS3FxsWY2m7VVq1ZpBw4c0G677TYtJiZGKy8v17tpIWPu3LladHS0Vlpaqp08edL9qqurc29zxx13aH379tX+8Y9/aJ9++qmWnZ2tZWdn69jq0NX0aRpN47X3lz179mgGg0F7/PHHtX/961/a66+/rtlsNm3NmjXubZ588kktJiZGe+edd7TPP/9cu+aaa/h4qQ/MmjVLS05Odj/au2HDBi0uLk6777773Nvw2gensA0jmqZpf/rTn7S+fftqJpNJGzVqlPbRRx/p3aSQAqDF18qVK93b1NfXa3feeafWs2dPzWazaddee6128uRJ/Rodws4PI7z2/vN///d/2uDBgzWz2awNHDhQ+8tf/uKxXlVVbeHChZrdbtfMZrN21VVXaYcOHdKptaGjurpau+eee7S+fftqFotFu/DCC7U//OEPmsPhcG/Dax+cBE1rMjUdERERUYCF5ZgRIiIiCh4MI0RERKQrhhEiIiLSFcMIERER6YphhIiIiHTFMEJERES6YhghIiIiXTGMEBERka4YRoioQ0pLSyEIQrPier62atUqxMTEuL9++OGHMXz4cL+ek4j0xTBCRC0aN24c5s+f7/569OjROHnyJKKjowPajgULFngUNiOi0GPQuwFE1D2YTCZdyqxHREQgIiIi4OclosBhzwgRNTN79mxs374dS5YsgSAIEAQBq1at8rhNc/Z2ysaNG5Geng6bzYbrrrsOdXV1WL16NVJTU9GzZ0/cfffdUBTFfWyHw4EFCxYgOTkZPXr0QFZWFkpLS1tty/m3aWbPno2pU6fi2WefRe/evREbG4u8vDy4XK5On4OI9MWeESJqZsmSJTh8+DAGDx6MRx99FADw1VdfNduurq4OL774IoqLi1FTU4Np06bh2muvRUxMDDZv3oxvvvkG06dPx5gxYzBjxgwAwLx583DgwAEUFxcjKSkJb731FiZOnIgvvvgC/fv371D73n//ffTu3Rvvv/8+jhw5ghkzZmD48OGYM2eOz85BRIHDMEJEzURHR8NkMsFms7lvzXz99dfNtnO5XFi2bBkuuugiAMB1112H1157DRUVFYiIiEBGRgbGjx+P999/HzNmzEBZWRlWrlyJsrIyJCUlAWgcE7J161asXLkSTzzxRIfa17NnT/z5z3+GJEkYOHAgJk+ejJKSEsyZM8dn5yCiwGEYIaJOs9ls7iACAHa7HampqR5jPOx2OyorKwEAX3zxBRRFwYABAzyO43A4EBsb2+HzXnzxxZAkyf1179698cUXX/j0HEQUOAwjRNRpRqPR42tBEFpcpqoqAODMmTOQJAl79+71CBMAvBqkGohzEFHgMIwQUYtMJpPHwFNfGDFiBBRFQWVlJa644gqfHjuQ5yAi3+LTNETUotTUVHz88cc4fvw4qqqq3D0PXTFgwADMnDkTN998MzZs2IBjx45hz549KCoqwqZNm3zQ6sCcg4h8i2GEiFq0YMECSJKEjIwMxMfHo6yszCfHXblyJW6++Wb87ne/Q3p6OqZOnYpPPvkEffv29cnxA3UOIvIdQdM0Te9GEBERUfhizwgRERHpimGEiIiIdMUwQkRERLpiGCEiIiJdMYwQERGRrhhGiIiISFcMI0RERKQrhhEiIiLSFcMIERER6YphhIiIiHTFMEJERES6+v81r419/70nKQAAAABJRU5ErkJggg==",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "kmf2.plot_survival_function(label=\"KM_estimate-2cm\")\n",
    "kmf5.plot_survival_function(label=\"KM_estimate-5cm\");"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### log-rank test\n",
    "\n",
    "- $ H_0 $ : the survival curves of the two groups are the same.\n",
    "- $ H_1 $ : the survival curves of the two groups are different."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>t_0</th>\n",
       "      <td>-1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>null_distribution</th>\n",
       "      <td>chi squared</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>degrees_of_freedom</th>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>test_name</th>\n",
       "      <td>logrank_test</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div><table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>test_statistic</th>\n",
       "      <th>p</th>\n",
       "      <th>-log2(p)</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>6.35</td>\n",
       "      <td>0.01</td>\n",
       "      <td>6.41</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>"
      ],
      "text/latex": [
       "\\begin{tabular}{lrrr}\n",
       " & test_statistic & p & -log2(p) \\\\\n",
       "0 & 6.35 & 0.01 & 6.41 \\\\\n",
       "\\end{tabular}\n"
      ],
      "text/plain": [
       "<lifelines.StatisticalResult: logrank_test>\n",
       "               t_0 = -1\n",
       " null_distribution = chi squared\n",
       "degrees_of_freedom = 1\n",
       "         test_name = logrank_test\n",
       "\n",
       "---\n",
       " test_statistic    p  -log2(p)\n",
       "           6.35 0.01      6.41"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "from lifelines.statistics import logrank_test\n",
    "\n",
    "res = logrank_test(T2, T5, E2, E5)\n",
    "\n",
    "res"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$ p < 0.05 $ ，按 $ \\alpha = 0.05 $ 水准，拒绝 $ H_0 $ ，接受 $ H_1 $ ，可以认为两组患者的术后生存率不同。"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": ".venv",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.10.12"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
